# -*- coding: utf-8 -*-
import numpy as np
import scipy.io as sio
import xlrd
import time
from matplotlib import pyplot as plt
import matplotlib
import h5py
from scipy.optimize import curve_fit
from scipy.ndimage import generic_filter

import os
import socket
os.system('chcp 65001')
Host_Name = socket.gethostname()
if Host_Name == 'VictorChengPC':
    Computer_Name = 'X:\\Work\\mail.ustc.edu.cn\\Rb87BEC - Files\\'
elif Host_Name == 'calcenter2':
    Computer_Name = 'D:\\vc\\mail.ustc.edu.cn\\Rb87BEC - Files\\'
font = {'size': 26}
matplotlib.rc('font', **font)
# %%
time_start = time.time()
time_start_struct = time.localtime(time_start)


# %%
def m_data(sp, **kw):
    if 'Find_Center' in kw and kw['Find_Center'] == 1:
        find_center_size = kw['Size']
        sp_mean = sp  # generic_filter(sp, np.nanmean, size=find_center_size, mode='mirror')
        return np.min(sp_mean), np.argmin(sp_mean), sp_mean
    size = kw['Size']
    # c1 = sp[0:size[0], 0:size[1]]
    # c2 = sp[-size[0]:, 0:size[1]]
    # c3 = sp[0:size[0], -size[1]:]
    # c4 = sp[-size[0]:, -size[1]:]
    c1 = sp[0, 0]
    c2 = sp[-1, 0]
    c3 = sp[0, -1]
    c4 = sp[-1, -1]
    return (np.mean(c1) + np.mean(c2) + np.mean(c3) + np.mean(c4)) / 4


def gamma_data(sp, **kw):
    if 'Find_Center' in kw and kw['Find_Center'] == 1:
        find_center_size = kw['Size']
        sp_mean = sp  # generic_filter(sp, np.nanmean, size=find_center_size, mode='mirror')
        return np.max(sp_mean), np.argmax(sp_mean), sp_mean
    data_size = np.shape(sp)
    size = kw['Size']
    # c = sp[data_size[0] // 2 - size[0]:data_size[0] // 2 + size[0],
    #        data_size[1] // 2 - size[1]:data_size[1] // 2 + size[1]]
    c = sp[data_size[0] // 2, data_size[1] // 2]
    return np.mean(c)


def x1(sp, **kw):
    if 'Find_Center' in kw and kw['Find_Center'] == 1:
        print("Not supported\n")
        return 0
    sp_mean = sp
    data_size = np.shape(sp)
    return sp_mean[data_size[0] // 2, 0]


def x2(sp, **kw):
    if 'Find_Center' in kw and kw['Find_Center'] == 1:
        print("Not supported\n")
        return 0
    sp_mean = sp
    data_size = np.shape(sp)
    return sp_mean[0, data_size[1] // 2]


# def gamma_data(sp, **kw):
#     sp_mean = sp
#     data_size = np.shape(sp)
#     return sp_mean[data_size[1] // 2, data_size[1] // 2]


# def m_data(sp, **kw):
#     sp_mean = sp
#     return sp_mean[0, 0]


def linear(x, a, b):
    return a * x + b


# %% Init
Fit_Prefix = ['NoFit', '']
Fit_BEC_Switch = 1
Remove_Upper_Switch = 1
Classify_Size = (8, 8)
Atom_Temperature = 150e-9
Data_Table = xlrd.open_workbook(Computer_Name + 'Rb87files\\3DSOC_spin_texture\\spin_texture.xlsx').sheets()[0]
Meta_Table = xlrd.open_workbook(Computer_Name + 'Rb87files\\3DSOC_spin_texture\\spin_texture.xlsx').sheets()[1]

# %% Calculate average
# @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
mz_simed_out = 0.25
if mz_simed_out == -0.25:
    Number_Remain = 24
    Meta_Row_List = range(202, 222 + 1)  # detuning = -0.5
elif mz_simed_out == 0:
    Number_Remain = 18
    Meta_Row_List = np.concatenate((np.arange(13, 51 + 1), np.arange(63, 66 + 1)))  # detuning = 0
elif mz_simed_out == -0.2:
    # Number_Remain = 18
    Number_Remain = 1000
    # Meta_Row_List = np.concatenate((np.arange(33, 71 + 1), np.arange(217, 217 + 1)))  # detuning = -0.4
    Meta_Row_List = np.arange(202, 222 + 1)
elif mz_simed_out == 0.25:
    Number_Remain = 1000
    Meta_Row_List = np.concatenate((np.arange(185, 192 + 1), np.arange(156, 168 + 1)))
# @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
Actual_Classify_All = []
Ave_Classify_All = []
Std_Classify_All = []
Ave_X1_All = []
Std_X1_All = []
Ave_X2_All = []
Std_X2_All = []
Ave_M_All = []
Std_M_All = []
Ave_Gamma_All = []
Std_Gamma_All = []
Ave_X1_Lower_All = []
Std_X1_Lower_All = []
Ave_X2_Lower_All = []
Std_X2_Lower_All = []
Ave_M_Lower_All = []
Std_M_Lower_All = []
Ave_Gamma_Lower_All = []
Std_Gamma_Lower_All = []

Detunings = []

Ave_Lower_Classify_All = []
Std_Lower_Classify_All = []
LenSPs = []
for Meta_Row_Index in range(len(Meta_Row_List)):
    Meta_Row = Meta_Row_List[Meta_Row_Index]
    Meta = Meta_Table.row_values(Meta_Row - 1)
    Try = int(Meta[1])
    Detuning = Meta[2]
    Detunings.append(Detuning)
    Mat_Data_command = Meta[3]
    Mat_Data = []
    exec(Mat_Data_command)
    Bad_Dict_Command = Meta[4]
    mz_simed = Meta[5]
    Bad_Dict = {}
    exec(Bad_Dict_Command)

    # Mat_Data = [593]
    # Bad_Dict = {}

    # Selected_Col = 83
    # Mat_Data = [364, 400, 414, 452, 478, 491, 539, 550]
    # Bad_Dict = {364:[1173], 539:[1260], 550:[1480], 478:[348], 400:[473]}

    Indicators = {}
    Indicators_Lower = {}
    SPs = {}
    SPs_Lower = {}
    for Line_Index in range(len(Mat_Data)):
        Line = Mat_Data[Line_Index]
        Input_Matrix = Data_Table.row_values(Line - 1)
        YearMonthDay = str(int(Input_Matrix[0]))
        Data_Folder_Name = Computer_Name + 'Rb87Data\\3DSOC_spin_texture\\' + \
            Fit_Prefix[Fit_BEC_Switch] + '\\' + YearMonthDay + '\\All_Spin_Polarization\\'
        Figure_Folder_Name = '\\\\calcenter2\\Rb87figures\\3DSOC_spin_texture\\figures' + \
                             Fit_Prefix[Fit_BEC_Switch] + '\\'+YearMonthDay + '\\All_Spin_Polarization\\'
        Ini_Photo = int(Input_Matrix[1])
        Fin_Photo = int(Input_Matrix[2])
        if Line in Bad_Dict.keys():
            Bad_List = Bad_Dict[Line]
        else:
            Bad_List = []
        for Photo_Count in range(Fin_Photo - Ini_Photo + 1):
            Photo_Ind = Photo_Count + Ini_Photo
            if Photo_Ind in Bad_List:
                continue
            Data = np.load(Data_Folder_Name + 'Cam_22860_photo%d_Spin_Polarization_FBZ_WithLower_%.3eK.npz' %
                           (Photo_Ind, Atom_Temperature))
            if Remove_Upper_Switch == 0:
                Data = np.load(Data_Folder_Name + 'Cam_22860_photo%d_Spin_Polarization_FBZ.npz' % Photo_Ind)
            Data_dict = dict(Data)
            SP = Data['SP']
            if 'SP_Lower' not in Data_dict:
                FBZ_Up_Lower = Data['FBZ_Up_Lower']
                FBZ_Dn_Lower = Data['FBZ_Dn_Lower']
                SP_Lower = (FBZ_Up_Lower - FBZ_Dn_Lower) / (FBZ_Up_Lower + FBZ_Dn_Lower)
                Data_dict['SP_Lower'] = SP_Lower
                np.savez(Data_Folder_Name + 'Cam_22860_photo%d_Spin_Polarization_FBZ_WithLower_%.3eK.npz' %
                         (Photo_Ind, Atom_Temperature), **Data_dict)
            else:
                SP_Lower = Data_dict['SP_Lower']
            if 'SP_Mean' not in Data_dict:
                SP_Mean = generic_filter(SP, np.nanmean, size=Classify_Size, mode='mirror')
                Data_dict['SP_Mean'] = SP_Mean
                np.savez(Data_Folder_Name + 'Cam_22860_photo%d_Spin_Polarization_FBZ_WithLower_%.3eK.npz' %
                         (Photo_Ind, Atom_Temperature), **Data_dict)
            else:
                SP_Mean = Data_dict['SP_Mean']
            if 'SP_Lower_Mean' not in Data_dict:
                SP_Lower_Mean = generic_filter(SP_Lower, np.nanmean, size=Classify_Size, mode='mirror')
                Data_dict['SP_Lower_Mean'] = SP_Lower_Mean
                np.savez(Data_Folder_Name + 'Cam_22860_photo%d_Spin_Polarization_FBZ_WithLower_%.3eK.npz' %
                         (Photo_Ind, Atom_Temperature), **Data_dict)
            else:
                SP_Lower_Mean = Data_dict['SP_Lower_Mean']

            M_Data = m_data(SP_Mean, Size=Classify_Size)
            Gamma_Data = gamma_data(SP_Mean, Size=Classify_Size)
            X1 = x1(SP_Mean)
            X2 = x2(SP_Mean)
            M_Data_Lower = m_data(SP_Lower_Mean, Size=Classify_Size)
            Gamma_Data_Lower = gamma_data(SP_Lower_Mean, Size=Classify_Size)
            X1_Lower = x1(SP_Lower_Mean)
            X2_Lower = x2(SP_Lower_Mean)

            if Detuning < 0:
                Indicator = M_Data
                Indicator_Lower = M_Data_Lower
                Classify = m_data
            else:
                Indicator = Gamma_Data
                Indicator_Lower = Gamma_Data_Lower
                Classify = gamma_data

            Indicators[(Line, Photo_Ind)] = Indicator
            Indicators_Lower[(Line, Photo_Ind)] = Indicator_Lower
            SPs[(Line, Photo_Ind)] = SP
            SPs_Lower[(Line, Photo_Ind)] = SP_Lower
            if len(SPs) >= Number_Remain:
                break
        if len(SPs) >= Number_Remain:
            break
    bins = np.linspace(-1, 1, 41)
    # fig = plt.figure(Meta_Row_Index)
    # hist_res = np.histogram(Indicators.values(), bins=bins)
    # plt.title('Meta_Row=%d, detuning=%.3f' % (Meta_Row, Detuning))
    # plt.show(block=False)
    # plt.pause(0.5)
    # plt.close(Meta_Row_Index)
    SP_Ave = np.average(np.array(list(SPs.values())), axis=0)
    SP_Std = np.std(np.array(list(SPs.values())), axis=0)

    np.savez(Computer_Name + "Rb87Data\\3DSOC_spin_texture\\Selected\\"
                             "ClassifyMeta_Row=%d,Detuning=%.3f.npz" % (Meta_Row, Detuning),
             Indicators=Indicators, Mat_Data=Mat_Data, Bad_Dict=Bad_Dict, SP_Ave=SP_Ave, SP_Std=SP_Std)
    sio.savemat(Computer_Name + "Rb87Data\\3DSOC_spin_texture\\Selected\\"
                                "ClassifyMeta_Row=%d,Detuning=%.3f.mat" % (Meta_Row, Detuning),
                {'Indicators': Indicators, 'Mat_Data': Mat_Data, 'Bad_Dict': Bad_Dict,
                    'SP_Ave': SP_Ave, 'SP_Std': SP_Std})

    SP_Ave_Mean = generic_filter(SP_Ave, np.nanmean, size=Classify_Size, mode='mirror')
    SP_Var_Mean = generic_filter(SP_Std ** 2, np.nanmean, size=Classify_Size, mode='mirror')

    Classify_Result = Classify(SP_Ave_Mean, Size=Classify_Size, Find_Center=1)
    Actual_Classify = np.unravel_index(Classify_Result[1], SP_Ave.shape)
    Actual_Classify_All.append(Actual_Classify)
    Ave_Classify_All.append(Classify_Result[0])
    Std_Classify_All.append(np.sqrt(SP_Var_Mean[Actual_Classify]))
    Data_Size = np.shape(SP_Std)

    X1_Result = x1(SP_Ave_Mean, Size=Classify_Size)
    Ave_X1_All.append(X1_Result)
    Std_X1_All.append(np.sqrt(SP_Var_Mean[Data_Size[0] // 2, 0]))
    X2_Result = x2(SP_Ave_Mean, Size=Classify_Size)
    Ave_X2_All.append(X2_Result)
    Std_X2_All.append(np.sqrt(SP_Var_Mean[0, Data_Size[1] // 2]))
    # M_Result = m_data(SP_Ave_Mean, Size=Classify_Size)
    if Detuning < 0:
        M_Result = SP_Ave_Mean[0, 0]
    else:
        M_Result = SP_Ave_Mean[0, 0]
    Ave_M_All.append(M_Result)
    Std_M_All.append(np.sqrt(SP_Var_Mean[0, 0]))
    # Gamma_Result = gamma_data(SP_Ave_Mean, Size=Classify_Size)
    Gamma_Result = SP_Ave_Mean[Data_Size[0] // 2, Data_Size[1] // 2]
    Ave_Gamma_All.append(Gamma_Result)
    Std_Gamma_All.append(np.sqrt(SP_Var_Mean[Data_Size[0] // 2, Data_Size[1] // 2]))

    SP_Lower_Ave = np.average(np.array(list(SPs_Lower.values())), axis=0)
    SP_Lower_Std = np.std(np.array(list(SPs_Lower.values())), axis=0)

    np.savez(Computer_Name + "Rb87Data\\3DSOC_spin_texture\\Selected\\"
                             "ClassifyMeta_Row=%d,Detuning=%.3f_Lower.npz" % (Meta_Row, Detuning),
             Indicators_Lower=Indicators_Lower, Mat_Data=Mat_Data, Bad_Dict=Bad_Dict,
             SP_Lower_Ave=SP_Lower_Ave, SP_Lower_Std=SP_Lower_Std)
    sio.savemat(Computer_Name + "Rb87Data\\3DSOC_spin_texture\\Selected\\"
                                "ClassifyMeta_Row=%d,Detuning=%.3f_Lower.mat" % (Meta_Row, Detuning),
                {'Indicators_Lower': Indicators_Lower, 'Mat_Data': Mat_Data, 'Bad_Dict': Bad_Dict,
                    'SP_Lower_Ave': SP_Lower_Ave, 'SP_Lower_Std': SP_Lower_Std})

    LenSPs.append(len(SPs))
    SP_Lower_Ave_Mean = generic_filter(SP_Lower_Ave, np.nanmean, size=Classify_Size, mode='mirror')
    SP_Lower_Var_Mean = generic_filter(SP_Lower_Std ** 2, np.nanmean, size=Classify_Size, mode='mirror')

    Classify_Lower_Result = Classify(SP_Lower_Ave_Mean, Size=Classify_Size, Find_Center=1)
    Actual_Classify_Lower = np.unravel_index(Classify_Lower_Result[1], SP_Lower_Ave.shape)
    Ave_Lower_Classify_All.append(Classify_Lower_Result[0])
    Std_Lower_Classify_All.append(np.sqrt(SP_Lower_Var_Mean[Actual_Classify_Lower]))

    X1_Lower_Result = x1(SP_Lower_Ave_Mean, Size=Classify_Size)
    Ave_X1_Lower_All.append(X1_Lower_Result)
    Std_X1_Lower_All.append(np.sqrt(SP_Lower_Var_Mean[Data_Size[0] // 2, 0]))
    X2_Lower_Result = x2(SP_Lower_Ave_Mean, Size=Classify_Size)
    Ave_X2_Lower_All.append(X2_Lower_Result)
    Std_X2_Lower_All.append(np.sqrt(SP_Lower_Var_Mean[0, Data_Size[1] // 2]))
    M_Lower_Result = m_data(SP_Lower_Ave_Mean, Size=Classify_Size)
    Ave_M_Lower_All.append(M_Lower_Result)
    Std_M_Lower_All.append(np.sqrt(SP_Lower_Var_Mean[Actual_Classify_Lower]))
    Gamma_Lower_Result = gamma_data(SP_Lower_Ave_Mean, Size=Classify_Size)
    Ave_Gamma_Lower_All.append(Gamma_Lower_Result)
    Std_Gamma_Lower_All.append(np.sqrt(SP_Lower_Var_Mean[Data_Size[0] // 2, Data_Size[1] // 2]))
    # fig.savefig("\\\\calcenter2\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\"
    #             "ClassifyMeta_Row=%d,Detuning=%.3f.png"%(Meta_Row, Detuning))
AveStd4All = np.zeros((len(Ave_M_All), 9))
AveStd4All[:, 0] = Ave_M_All
AveStd4All[:, 1] = Std_M_All
AveStd4All[:, 2] = Ave_X1_All
AveStd4All[:, 3] = Std_X1_All
AveStd4All[:, 4] = Ave_X2_All
AveStd4All[:, 5] = Std_X2_All
AveStd4All[:, 6] = Ave_Gamma_All
AveStd4All[:, 7] = Std_Gamma_All

AveStd4All_Lower = np.zeros((len(Ave_M_Lower_All), 9))
AveStd4All_Lower[:, 0] = Ave_M_Lower_All
AveStd4All_Lower[:, 1] = Std_M_Lower_All
AveStd4All_Lower[:, 2] = Ave_X1_Lower_All
AveStd4All_Lower[:, 3] = Std_X1_Lower_All
AveStd4All_Lower[:, 4] = Ave_X2_Lower_All
AveStd4All_Lower[:, 5] = Std_X2_Lower_All
AveStd4All_Lower[:, 6] = Ave_Gamma_Lower_All
AveStd4All_Lower[:, 7] = Std_Gamma_Lower_All

Data_Labels = ['M_Ave', 'M_Std', 'X1_Ave', 'X1_Std', 'X2_Ave', 'X2_Std', 'Gamma_Ave', 'Gamma_Std', 'kz']
AveStd4All_Save = AveStd4All.copy()
AveStd4All_Lower_Save = AveStd4All_Lower.copy()
#%% -0.4Er
AveStd4All = AveStd4All_Save.copy()
AveStd4All_Lower = AveStd4All_Lower_Save.copy()
if mz_simed_out == -0.2 or mz_simed_out == 0 or mz_simed_out == 0.25:
    Pair_Data = sio.loadmat(Computer_Name + "Rb87Data\\3DSOCbandData\\200117\\"
                                            "detuning_ky_pair,mz_sim=%.3f.mat" % mz_simed_out)
    Detuning_ky_pair = Pair_Data['detuning_ky_pair']
    if mz_simed_out == -0.2:
        Detuning_ky_pair[0, 1] = 1
        Detuning_ky_pair[-1, 1] = 0
    for_del = np.asarray(Detuning_ky_pair[:, 1] == 100).nonzero()
    Detuning_ky_pair = np.delete(Detuning_ky_pair, for_del, 0)
    Detunings = Detuning_ky_pair[:, 0]
    kys = Detuning_ky_pair[:, 1]
    AveStd4All[:, 8] = kys
    kys_sort_arg = list(np.argsort(-kys))

    LenSPs = np.array(LenSPs)
    LenSPs = LenSPs[kys_sort_arg]
    kys = kys[kys_sort_arg]
    Detunings = Detunings[kys_sort_arg]
    AveStd4All = AveStd4All[kys_sort_arg, :]

    Ave_M_All = AveStd4All[:, 0]
    Std_M_All = AveStd4All[:, 1]
    Ave_X1_All = AveStd4All[:, 2]
    Std_X1_All = AveStd4All[:, 3]
    Ave_X2_All = AveStd4All[:, 4]
    Std_X2_All = AveStd4All[:, 5]
    Ave_Gamma_All = AveStd4All[:, 6]
    Std_Gamma_All = AveStd4All[:, 7]
    kys = AveStd4All[:, 8]

    ky_flipped = np.concatenate((kys, -np.flip(kys)))
    Ave_M_All_flipped = np.concatenate((Ave_M_All, np.flip(Ave_M_All)))
    Std_M_All_flipped = np.concatenate((Std_M_All, np.flip(Std_M_All)))
    Ave_X1_All_flipped = np.concatenate((Ave_X1_All, np.flip(Ave_X1_All)))
    Std_X1_All_flipped = np.concatenate((Std_X1_All, np.flip(Std_X1_All)))
    Ave_X2_All_flipped = np.concatenate((Ave_X2_All, np.flip(Ave_X2_All)))
    Std_X2_All_flipped = np.concatenate((Std_X2_All, np.flip(Std_X2_All)))
    Ave_Gamma_All_flipped = np.concatenate((Ave_Gamma_All, np.flip(Ave_Gamma_All)))
    Std_Gamma_All_flipped = np.concatenate((Std_Gamma_All, np.flip(Std_Gamma_All)))

    fig = plt.figure(figsize=(20, 15))
    # plt.plot(ky_flipped, Ave_Classify_All_Fit_flipped, 'r-', linewidth=5)
    plt.errorbar(ky_flipped, Ave_M_All_flipped, Std_M_All_flipped,
                 fmt='bo-', label='M', linewidth=5, markersize=20, mec='k', mfc='k')
    plt.errorbar(ky_flipped, Ave_X1_All_flipped, Std_X1_All_flipped,
                 fmt='mo-', label=r'$X_1$', linewidth=5, markersize=20, mec='k', mfc='k')
    plt.errorbar(ky_flipped, Ave_X2_All_flipped, Std_X2_All_flipped,
                 fmt='yo-', label=r'$X_2$', linewidth=5, markersize=20, mec='k', mfc='k')
    plt.errorbar(ky_flipped, Ave_Gamma_All_flipped, Std_Gamma_All_flipped,
                 fmt='go-', label=r'$\Gamma$', linewidth=5, markersize=20, mec='k', mfc='k')

    # plt.plot((np.min(ky_flipped), np.max(ky_flipped)), (0, 0))
    # plt.title('Average Photo Center Value')
    plt.xlabel(r'kz/$\pi$')
    plt.ylabel(r'$P(\Lambda_i)$')
    plt.xlim([np.min(ky_flipped), np.max(ky_flipped)])
    plt.ylim([-1, 1])

    # plt.plot(detuning_PWM_Weyl, 0, 'o', ms=10)
    # plt.title('Ave Photo Cntr Value with upper, Weyl:detuning=%.3f, T=%.3eK(SD)' %
    #           (detuning_PWM_Weyl, Atom_Temperature))
    plt.legend(loc='upper right')
    plt.tick_params(labelsize=30)
    plt.show()
    if Host_Name == 'calcenter2':
        fig.savefig('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value with upper, mz_simed_out=%.3f, T=%.3eK(SD).png' %
                    (mz_simed_out, Atom_Temperature))
    elif Host_Name == 'VictorChengPC':
        fig.savefig('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value with upper, mz_simed_out=%.3f, T=%.3eK(SD).png' %
                    (mz_simed_out, Atom_Temperature))

    LenSPs = np.array(LenSPs)
    LenSPs_flipped = np.concatenate((LenSPs, np.flip(LenSPs)))
    fig = plt.figure(figsize=(20, 15))
    # plt.plot(ky_flipped, Ave_Classify_All_Fit_flipped, 'r-', linewidth=5)
    plt.errorbar(ky_flipped, Ave_M_All_flipped, Std_M_All_flipped / np.sqrt(LenSPs_flipped),
                 fmt='bo-', label='M', linewidth=5, markersize=15, mec='k', mfc='k')
    plt.errorbar(ky_flipped, Ave_X1_All_flipped, Std_X1_All_flipped / np.sqrt(LenSPs_flipped),
                 fmt='mo-', label=r'$X_1$', linewidth=5, markersize=15, mec='k', mfc='k')
    plt.errorbar(ky_flipped, Ave_X2_All_flipped, Std_X2_All_flipped / np.sqrt(LenSPs_flipped),
                 fmt='yo-', label=r'$X_2$', linewidth=5, markersize=15, mec='k', mfc='k')
    plt.errorbar(ky_flipped, Ave_Gamma_All_flipped, Std_Gamma_All_flipped / np.sqrt(LenSPs_flipped),
                 fmt='go-', label=r'$\Gamma$', linewidth=5, markersize=15, mec='k', mfc='k')

    plt.plot((np.min(ky_flipped), np.max(ky_flipped)), (0, 0))
    # plt.title('Average Photo Center Value')
    plt.xlabel(r'kz/$\pi$')
    plt.ylabel(r'$P(\Lambda_i)$')
    plt.xlim([np.min(ky_flipped), np.max(ky_flipped)])
    plt.ylim([-1, 1])
    # plt.plot(detuning_PWM_Weyl, 0, 'o', ms=10)
    # plt.title('Ave Photo Cntr Value with upper, Weyl:detuning=%.3f, T=%.3eK(SE)' %
    #           (detuning_PWM_Weyl, Atom_Temperature))
    plt.legend(loc='upper right')
    plt.tick_params(labelsize=30)
    plt.show()
    if Host_Name == 'calcenter2':
        fig.savefig('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value with upper, mz_simed_out=%.3f, T=%.3eK(SE).png' %
                    (mz_simed_out, Atom_Temperature))
    elif Host_Name == 'VictorChengPC':
        fig.savefig('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value with upper, mz_simed_out=%.3f, T=%.3eK(SE).png' %
                    (mz_simed_out, Atom_Temperature))
    if Host_Name == 'calcenter2':
        np.savez('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                 'Ave Photo Cntr Value with upper, mz_simed_out=%.3f, T=%.3eKData.npz' %
                 (mz_simed_out, Atom_Temperature), AveStd4All=AveStd4All, Data_Labels=Data_Labels)
        sio.savemat('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value with upper, mz_simed_out=%.3f, T=%.3eKData.mat' %
                    (mz_simed_out, Atom_Temperature), {'AveStd4All': AveStd4All, 'Data_Labels': Data_Labels})
    elif Host_Name == 'VictorChengPC':
        np.savez('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                 'Ave Photo Cntr Value with upper, mz_simed_out=%.3f, T=%.3eKData.npz' %
                 (mz_simed_out, Atom_Temperature), AveStd4All=AveStd4All, Data_Labels=Data_Labels)
        sio.savemat('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value with upper, mz_simed_out=%.3f, T=%.3eKData.mat' %
                    (mz_simed_out, Atom_Temperature), {'AveStd4All': AveStd4All, 'Data_Labels': Data_Labels})

    # Lower

    # Pair_Data = sio.loadmat(Computer_Name + "Rb87Data\\3DSOCbandData\\190827\\"
    #                                         "detuning_ky_pair,mz_sim=%.3f.mat" % mz_simed_out)
    Detuning_ky_pair = Pair_Data['detuning_ky_pair']
    if mz_simed_out == -0.2:
        Detuning_ky_pair[0, 1] = 1
        Detuning_ky_pair[-1, 1] = 0
    for_del = np.asarray(Detuning_ky_pair[:, 1] == 100).nonzero()
    Detuning_ky_pair = np.delete(Detuning_ky_pair, for_del, 0)
    Detunings = Detuning_ky_pair[:, 0]
    kys = Detuning_ky_pair[:, 1]
    AveStd4All_Lower[:, 8] = kys
    kys_sort_arg = list(np.argsort(-kys))

    kys = kys[kys_sort_arg]
    Detunings = Detunings[kys_sort_arg]
    AveStd4All_Lower = AveStd4All_Lower[kys_sort_arg, :]

    Ave_M_Lower_All = AveStd4All_Lower[:, 0]
    Std_M_Lower_All = AveStd4All_Lower[:, 1]
    Ave_X1_Lower_All = AveStd4All_Lower[:, 2]
    Std_X1_Lower_All = AveStd4All_Lower[:, 3]
    Ave_X2_Lower_All = AveStd4All_Lower[:, 4]
    Std_X2_Lower_All = AveStd4All_Lower[:, 5]
    Ave_Gamma_Lower_All = AveStd4All_Lower[:, 6]
    Std_Gamma_Lower_All = AveStd4All_Lower[:, 7]
    kys = AveStd4All_Lower[:, 8]

    ky_flipped = np.concatenate((kys, -np.flip(kys)))
    Ave_M_Lower_All_flipped = np.concatenate((Ave_M_Lower_All, np.flip(Ave_M_Lower_All)))
    Std_M_Lower_All_flipped = np.concatenate((Std_M_Lower_All, np.flip(Std_M_Lower_All)))
    Ave_X1_Lower_All_flipped = np.concatenate((Ave_X1_Lower_All, np.flip(Ave_X1_Lower_All)))
    Std_X1_Lower_All_flipped = np.concatenate((Std_X1_Lower_All, np.flip(Std_X1_Lower_All)))
    Ave_X2_Lower_All_flipped = np.concatenate((Ave_X2_Lower_All, np.flip(Ave_X2_Lower_All)))
    Std_X2_Lower_All_flipped = np.concatenate((Std_X2_Lower_All, np.flip(Std_X2_Lower_All)))
    Ave_Gamma_Lower_All_flipped = np.concatenate((Ave_Gamma_Lower_All, np.flip(Ave_Gamma_Lower_All)))
    Std_Gamma_Lower_All_flipped = np.concatenate((Std_Gamma_Lower_All, np.flip(Std_Gamma_Lower_All)))

    fig = plt.figure(figsize=(20, 15))
    plt.errorbar(ky_flipped, Ave_M_Lower_All_flipped, Std_M_Lower_All_flipped,
                 fmt='bo-', label='M', linewidth=5, markersize=20, mec='k', mfc='k')
    plt.errorbar(ky_flipped, Ave_X1_Lower_All_flipped, Std_X1_Lower_All_flipped,
                 fmt='mo-', label=r'$X_1$', linewidth=5, markersize=20, mec='k', mfc='k')
    plt.errorbar(ky_flipped, Ave_X2_Lower_All_flipped, Std_X2_Lower_All_flipped,
                 fmt='yo-', label=r'$X_2$', linewidth=5, markersize=20, mec='k', mfc='k')
    plt.errorbar(ky_flipped, Ave_Gamma_Lower_All_flipped, Std_Gamma_Lower_All_flipped,
                 fmt='go-', label=r'$\Gamma$', linewidth=5, markersize=20, mec='k', mfc='k')
    plt.plot((np.min(ky_flipped), np.max(ky_flipped)), (0, 0))
    plt.xlabel(r'kz/$\pi$')
    plt.ylabel(r'$P(\Lambda_i)$')
    plt.xlim([np.min(ky_flipped), np.max(ky_flipped)])
    plt.ylim([-1, 1])
    plt.legend(loc='upper right')
    plt.tick_params(labelsize=30)
    plt.show()

    fig = plt.figure(figsize=(20, 15))
    plt.errorbar(ky_flipped, Ave_M_Lower_All_flipped, Std_M_Lower_All_flipped / np.sqrt(LenSPs_flipped),
                 fmt='bo-', label='M', linewidth=5, markersize=15, mec='k', mfc='k')
    plt.errorbar(ky_flipped, Ave_X1_Lower_All_flipped, Std_X1_Lower_All_flipped / np.sqrt(LenSPs_flipped),
                 fmt='mo-', label=r'$X_1$', linewidth=5, markersize=15, mec='k', mfc='k')
    plt.errorbar(ky_flipped, Ave_X2_Lower_All_flipped, Std_X2_Lower_All_flipped / np.sqrt(LenSPs_flipped),
                 fmt='yo-', label=r'$X_2$', linewidth=5, markersize=15, mec='k', mfc='k')
    plt.errorbar(ky_flipped, Ave_Gamma_Lower_All_flipped, Std_Gamma_Lower_All_flipped / np.sqrt(LenSPs_flipped),
                 fmt='go-', label=r'$\Gamma$', linewidth=5, markersize=15, mec='k', mfc='k')
    plt.plot((np.min(ky_flipped), np.max(ky_flipped)), (0, 0))
    plt.xlabel(r'kz/$\pi$')
    plt.ylabel(r'$P(\Lambda_i)$')
    plt.xlim([np.min(ky_flipped), np.max(ky_flipped)])
    plt.ylim([-1, 1])
    plt.legend(loc='upper right')
    plt.tick_params(labelsize=30)
    plt.show()
    if Host_Name == 'calcenter2':
        fig.savefig('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value no upper, mz_simed_out=%.3f, T=%.3eK(SE).png' %
                    (mz_simed_out, Atom_Temperature))
    elif Host_Name == 'VictorChengPC':
        fig.savefig('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value no upper, mz_simed_out=%.3f, T=%.3eK(SE).png' %
                    (mz_simed_out, Atom_Temperature))
    if Host_Name == 'calcenter2':
        np.savez('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                 'Ave Photo Cntr Value no upper, mz_simed_out=%.3f, T=%.3eKData.npz' %
                 (mz_simed_out, Atom_Temperature), AveStd4All_Lower=AveStd4All_Lower, Data_Labels=Data_Labels)
        sio.savemat('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value no upper, mz_simed_out=%.3f, T=%.3eKData.mat' %
                    (mz_simed_out, Atom_Temperature),
                    {'AveStd4All_Lower': AveStd4All_Lower, 'Data_Labels': Data_Labels})
    elif Host_Name == 'VictorChengPC':
        np.savez('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                 'Ave Photo Cntr Value no upper, mz_simed_out=%.3f, T=%.3eKData.npz' %
                 (mz_simed_out, Atom_Temperature), AveStd4All_Lower=AveStd4All_Lower, Data_Labels=Data_Labels)
        sio.savemat('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value no upper, mz_simed_out=%.3f, T=%.3eKData.mat' %
                    (mz_simed_out, Atom_Temperature),
                    {'AveStd4All_Lower': AveStd4All_Lower, 'Data_Labels': Data_Labels})
    # Select

    p_res, pcov = curve_fit(linear, Detunings, Ave_M_All, sigma=Std_M_All)
    Ave_M_All_Fits = p_res[0] * Detunings + p_res[1]

    mz_simed_out = -0.2
    if mz_simed_out == -0.25:
        Number_Remain = 24
        Meta_Row_List = range(202, 222 + 1)  # detuning = -0.5
    elif mz_simed_out == 0:
        Number_Remain = 18
        Meta_Row_List = np.concatenate((np.arange(13, 51 + 1), np.arange(63, 66 + 1)))  # detuning = 0
    elif mz_simed_out == -0.2:
        Number_Remain = 18
        Meta_Row_List = np.concatenate((np.arange(33, 71 + 1), np.arange(217, 217 + 1)))  # detuning = -0.4

    Actual_Classify_All_Sel = []
    Ave_Classify_All_Sel = []
    Std_Classify_All_Sel = []
    Ave_Lower_Classify_All_Sel = []
    Std_Lower_Classify_All_Sel = []
    AveStd4All_Sel = np.zeros((len(Meta_Row_List), 10))
    AveStd4All_Lower_Sel = np.zeros((len(Meta_Row_List), 10))
    LenSPs_Sel = []
    for Meta_Row_Index in range(len(Meta_Row_List)):
        Ave_M_All_Fit = Ave_M_All_Fits[Meta_Row_Index]
        Meta_Row = Meta_Row_List[Meta_Row_Index]
        Meta = Meta_Table.row_values(Meta_Row - 1)
        Try = int(Meta[1])
        Detuning = Meta[2]
        Mat_Data_command = Meta[3]
        Mat_Data = []
        exec(Mat_Data_command)
        Bad_Dict_Command = Meta[4]
        mz_simed = Meta[5]
        Bad_Dict = {}
        exec(Bad_Dict_Command)

        # Mat_Data = [593]
        # Bad_Dict = {}

        # Selected_Col = 83
        # Mat_Data = [364, 400, 414, 452, 478, 491, 539, 550]
        # Bad_Dict = {364:[1173], 539:[1260], 550:[1480], 478:[348], 400:[473]}

        Indicators_Sel = {}
        Indicators_Lower_Sel = {}
        SPs_Sel = {}
        SPs_Lower_Sel = {}
        for Line_Index in range(len(Mat_Data)):
            Line = Mat_Data[Line_Index]
            Input_Matrix = Data_Table.row_values(Line - 1)
            YearMonthDay = str(int(Input_Matrix[0]))
            Data_Folder_Name = Computer_Name + 'Rb87Data\\3DSOC_spin_texture\\' + \
                Fit_Prefix[Fit_BEC_Switch] + '\\' + YearMonthDay + '\\All_Spin_Polarization\\'
            Figure_Folder_Name = '\\\\calcenter2\\Rb87figures\\3DSOC_spin_texture\\figures' + \
                                 Fit_Prefix[Fit_BEC_Switch] + '\\'+YearMonthDay + '\\All_Spin_Polarization\\'
            Ini_Photo = int(Input_Matrix[1])
            Fin_Photo = int(Input_Matrix[2])
            if Line in Bad_Dict.keys():
                Bad_List = Bad_Dict[Line]
            else:
                Bad_List = []
            for Photo_Count in range(Fin_Photo - Ini_Photo + 1):
                Photo_Ind = Photo_Count + Ini_Photo
                if Photo_Ind in Bad_List:
                    continue
                Data = np.load(Data_Folder_Name + 'Cam_22860_photo%d_Spin_Polarization_FBZ_WithLower_%.3eK.npz' %
                               (Photo_Ind, Atom_Temperature))
                if Remove_Upper_Switch == 0:
                    Data = np.load(Data_Folder_Name + 'Cam_22860_photo%d_Spin_Polarization_FBZ.npz' % Photo_Ind)
                Data_dict = dict(Data)
                SP = Data['SP']
                if 'SP_Lower' not in Data_dict:
                    FBZ_Up_Lower = Data['FBZ_Up_Lower']
                    FBZ_Dn_Lower = Data['FBZ_Dn_Lower']
                    SP_Lower = (FBZ_Up_Lower - FBZ_Dn_Lower) / (FBZ_Up_Lower + FBZ_Dn_Lower)
                    Data_dict['SP_Lower'] = SP_Lower
                    np.savez(Data_Folder_Name + 'Cam_22860_photo%d_Spin_Polarization_FBZ_WithLower_%.3eK.npz' %
                             (Photo_Ind, Atom_Temperature), **Data_dict)
                else:
                    SP_Lower = Data_dict['SP_Lower']
                if 'SP_Mean' not in Data_dict:
                    SP_Mean = generic_filter(SP, np.nanmean, size=Classify_Size, mode='mirror')
                    Data_dict['SP_Mean'] = SP_Mean
                    np.savez(Data_Folder_Name + 'Cam_22860_photo%d_Spin_Polarization_FBZ_WithLower_%.3eK.npz' %
                             (Photo_Ind, Atom_Temperature), **Data_dict)
                else:
                    SP_Mean = Data_dict['SP_Mean']
                if 'SP_Lower_Mean' not in Data_dict:
                    SP_Lower_Mean = generic_filter(SP_Lower, np.nanmean, size=Classify_Size, mode='mirror')
                    Data_dict['SP_Lower_Mean'] = SP_Lower_Mean
                    np.savez(Data_Folder_Name + 'Cam_22860_photo%d_Spin_Polarization_FBZ_WithLower_%.3eK.npz' %
                             (Photo_Ind, Atom_Temperature), **Data_dict)
                else:
                    SP_Lower_Mean = Data_dict['SP_Lower_Mean']

                M_Data = m_data(SP_Mean, Size=Classify_Size)
                if np.abs(M_Data - Ave_M_All_Fit) >= 0.05:
                    continue  # @@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
                Gamma_Data = gamma_data(SP_Mean, Size=Classify_Size)
                X1 = x1(SP_Mean)
                X2 = x2(SP_Mean)
                M_Data_Lower = m_data(SP_Lower_Mean, Size=Classify_Size)
                Gamma_Data_Lower = gamma_data(SP_Lower_Mean, Size=Classify_Size)
                X1_Lower = x1(SP_Lower_Mean)
                X2_Lower = x2(SP_Lower_Mean)

                if Detuning < 0:
                    Indicator = M_Data
                    Indicator_Lower = M_Data_Lower
                    Classify = m_data
                else:
                    Indicator = Gamma_Data
                    Indicator_Lower = Gamma_Data_Lower
                    Classify = gamma_data

                Indicators_Sel[(Line, Photo_Ind)] = Indicator
                Indicators_Lower_Sel[(Line, Photo_Ind)] = Indicator_Lower
                SPs_Sel[(Line, Photo_Ind)] = SP
                SPs_Lower_Sel[(Line, Photo_Ind)] = SP_Lower
                if len(SPs_Sel) >= Number_Remain:
                    break
            if len(SPs_Sel) >= Number_Remain:
                break
        bins = np.linspace(-1, 1, 41)
        # fig = plt.figure(Meta_Row_Index)
        # hist_res = np.histogram(Indicators.values(), bins=bins)
        # plt.title('Meta_Row=%d, detuning=%.3f' % (Meta_Row, Detuning))
        # plt.show(block=False)
        # plt.pause(0.5)
        # plt.close(Meta_Row_Index)
        if len(SPs_Sel) == 1:
            # print('aha')
            pass
        if len(SPs_Sel) == 0:
            # print('oh shit')
            time.sleep(0.01)
            SPs_Sel[(Line, Photo_Ind)] = np.zeros_like(SP)
            SPs_Lower_Sel[(Line, Photo_Ind)] = np.zeros_like(SP)
            AveStd4All_Sel[Meta_Row_Index, 9] = 1
            AveStd4All_Lower_Sel[Meta_Row_Index, 9] = 1

        SP_Ave_Sel = np.average(np.array(list(SPs_Sel.values())), axis=0)
        SP_Std_Sel = np.std(np.array(list(SPs_Sel.values())), axis=0)

        np.savez(Computer_Name + "Rb87Data\\3DSOC_spin_texture\\Selected\\"
                                 "ClassifyMeta_Row=%d,Detuning=%.3f_Sel.npz" % (Meta_Row, Detuning),
                 Indicators_Sel=Indicators_Sel, Mat_Data=Mat_Data, Bad_Dict=Bad_Dict,
                 SP_Ave_Sel=SP_Ave_Sel, SP_Std_Sel=SP_Std_Sel)
        sio.savemat(Computer_Name + "Rb87Data\\3DSOC_spin_texture\\Selected\\"
                                    "ClassifyMeta_Row=%d,Detuning=%.3f_Sel.mat" % (Meta_Row, Detuning),
                    {'Indicators_Sel': Indicators_Sel, 'Mat_Data': Mat_Data, 'Bad_Dict': Bad_Dict,
                        'SP_Ave_Sel': SP_Ave_Sel, 'SP_Std_Sel': SP_Std_Sel})
        SP_Ave_Mean_Sel = generic_filter(SP_Ave_Sel, np.nanmean, size=Classify_Size, mode='mirror')
        SP_Var_Mean_Sel = generic_filter(SP_Std_Sel ** 2, np.nanmean, size=Classify_Size, mode='mirror')

        Classify_Result_Sel = Classify(SP_Ave_Mean_Sel, Size=Classify_Size, Find_Center=1)
        Actual_Classify_Sel = np.unravel_index(Classify_Result_Sel[1], SP_Ave_Sel.shape)
        Actual_Classify_All_Sel.append(Actual_Classify_Sel)
        Ave_Classify_All_Sel.append(Classify_Result_Sel[0])
        Std_Classify_All_Sel.append(np.sqrt(SP_Var_Mean_Sel[Actual_Classify_Sel]))
        Data_Size = np.shape(SP_Std_Sel)

        X1_Result_Sel = x1(SP_Ave_Mean_Sel, Size=Classify_Size)
        AveStd4All_Sel[Meta_Row_Index, 2] = X1_Result_Sel
        AveStd4All_Sel[Meta_Row_Index, 3] = np.sqrt(SP_Var_Mean_Sel[Data_Size[0] // 2, 0])

        X2_Result_Sel = x2(SP_Ave_Mean_Sel, Size=Classify_Size)
        AveStd4All_Sel[Meta_Row_Index, 4] = X2_Result_Sel
        AveStd4All_Sel[Meta_Row_Index, 5] = np.sqrt(SP_Var_Mean_Sel[0, Data_Size[1] // 2])
        # M_Result = m_data(SP_Ave_Mean, Size=Classify_Size)
        if Detuning < 0:
            M_Result_Sel = SP_Ave_Mean_Sel[0, 0]
        else:
            M_Result_Sel = SP_Ave_Mean_Sel[0, 0]
        AveStd4All_Sel[Meta_Row_Index, 0] = M_Result_Sel
        AveStd4All_Sel[Meta_Row_Index, 1] = np.sqrt(SP_Var_Mean_Sel[0, 0])
        Gamma_Result_Sel = SP_Ave_Mean_Sel[Data_Size[0] // 2, Data_Size[1] // 2]
        AveStd4All_Sel[Meta_Row_Index, 6] = Gamma_Result_Sel
        AveStd4All_Sel[Meta_Row_Index, 7] = np.sqrt(SP_Var_Mean_Sel[Data_Size[0] // 2, Data_Size[1] // 2])

        SP_Lower_Ave_Sel = np.average(np.array(list(SPs_Lower_Sel.values())), axis=0)
        SP_Lower_Std_Sel = np.std(np.array(list(SPs_Lower_Sel.values())), axis=0)

        np.savez(Computer_Name + "Rb87Data\\3DSOC_spin_texture\\Selected\\"
                                 "ClassifyMeta_Row=%d,Detuning=%.3f_Lower_Sel.npz" % (Meta_Row, Detuning),
                 Indicators_Lower_Sel=Indicators_Lower_Sel, Mat_Data=Mat_Data, Bad_Dict=Bad_Dict,
                 SP_Lower_Ave_Sel=SP_Lower_Ave_Sel, SP_Lower_Std_Sel=SP_Lower_Std_Sel)
        sio.savemat(Computer_Name + "Rb87Data\\3DSOC_spin_texture\\Selected\\"
                                    "ClassifyMeta_Row=%d,Detuning=%.3f_Lower.mat" % (Meta_Row, Detuning),
                    {'Indicators_Lower_Sel': Indicators_Lower_Sel, 'Mat_Data': Mat_Data, 'Bad_Dict': Bad_Dict,
                        'SP_Lower_Ave_Sel': SP_Lower_Ave_Sel, 'SP_Lower_Std_Sel': SP_Lower_Std_Sel})

        LenSPs_Sel.append(len(SPs_Sel))
        SP_Lower_Ave_Mean_Sel = generic_filter(SP_Lower_Ave_Sel, np.nanmean, size=Classify_Size, mode='mirror')
        SP_Lower_Var_Mean_Sel = generic_filter(SP_Lower_Std_Sel ** 2, np.nanmean, size=Classify_Size, mode='mirror')

        Classify_Lower_Result_Sel = Classify(SP_Lower_Ave_Mean_Sel, Size=Classify_Size, Find_Center=1)
        Actual_Classify_Lower_Sel = np.unravel_index(Classify_Lower_Result_Sel[1], SP_Lower_Ave_Sel.shape)
        Ave_Lower_Classify_All_Sel.append(Classify_Lower_Result_Sel[0])
        Std_Lower_Classify_All_Sel.append(np.sqrt(SP_Lower_Var_Mean_Sel[Actual_Classify_Lower_Sel]))

        X1_Lower_Result_Sel = x1(SP_Lower_Ave_Mean_Sel, Size=Classify_Size)
        AveStd4All_Lower_Sel[Meta_Row_Index, 2] = X1_Lower_Result_Sel
        AveStd4All_Lower_Sel[Meta_Row_Index, 3] = np.sqrt(SP_Lower_Var_Mean_Sel[Data_Size[0] // 2, 0])

        X2_Lower_Result_Sel = x2(SP_Lower_Ave_Mean_Sel, Size=Classify_Size)
        AveStd4All_Lower_Sel[Meta_Row_Index, 4] = X2_Lower_Result_Sel
        AveStd4All_Lower_Sel[Meta_Row_Index, 5] = np.sqrt(SP_Lower_Var_Mean_Sel[0, Data_Size[1] // 2])

        M_Lower_Result_Sel = m_data(SP_Lower_Ave_Mean_Sel, Size=Classify_Size)
        AveStd4All_Lower_Sel[Meta_Row_Index, 0] = M_Lower_Result_Sel
        AveStd4All_Lower_Sel[Meta_Row_Index, 1] = np.sqrt(SP_Lower_Var_Mean_Sel[0, 0])

        Gamma_Lower_Result_Sel = gamma_data(SP_Lower_Ave_Mean_Sel, Size=Classify_Size)
        AveStd4All_Lower_Sel[Meta_Row_Index, 6] = Gamma_Lower_Result_Sel
        AveStd4All_Lower_Sel[Meta_Row_Index, 7] = np.sqrt(SP_Lower_Var_Mean_Sel[Data_Size[0] // 2, Data_Size[1] // 2])

    # Seled, plot
    # Pair_Data = sio.loadmat(Computer_Name + "Rb87Data\\3DSOCbandData\\190827\\"
    #                                         "detuning_ky_pair,mz_sim=%.3f.mat" % mz_simed_out)
    Detuning_ky_pair = Pair_Data['detuning_ky_pair']
    if mz_simed_out == -0.2:
        Detuning_ky_pair[0, 1] = 1
        Detuning_ky_pair[-1, 1] = 0
    for_del = np.asarray(Detuning_ky_pair[:, 1] == 100).nonzero()
    Detuning_ky_pair = np.delete(Detuning_ky_pair, for_del, 0)
    Detunings = Detuning_ky_pair[:, 0]
    kys = Detuning_ky_pair[:, 1]
    AveStd4All_Sel[:, 8] = kys
    kys_sort_arg = list(np.argsort(-kys))
    LenSPs_Sel = np.array(LenSPs_Sel)

    LenSPs_Sel = LenSPs_Sel[kys_sort_arg]
    kys = kys[kys_sort_arg]
    Detunings = Detunings[kys_sort_arg]
    AveStd4All_Sel = AveStd4All_Sel[kys_sort_arg, :]

    ky_flipped = np.concatenate((kys, -np.flip(kys)))
    AveStd4All_Sel_flipped = np.concatenate((AveStd4All_Sel, np.flip(AveStd4All_Sel, axis=0)), axis=0)
    AveStd4All_Sel_flipped[:, 8] = ky_flipped
    for_del = np.asarray(AveStd4All_Sel_flipped[:, -1]).nonzero()
    if Host_Name == 'calcenter2':
        np.savez('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                 'Ave Photo Cntr Value with upper, mz_simed_out=%.3f, T=%.3eKData_Sel.npz' %
                 (mz_simed_out, Atom_Temperature), AveStd4All_Sel=AveStd4All_Sel, Data_Labels=Data_Labels)
        sio.savemat('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value with upper, mz_simed_out=%.3f, T=%.3eKData_Sel.mat' %
                    (mz_simed_out, Atom_Temperature), {'AveStd4All_Sel': AveStd4All_Sel, 'Data_Labels': Data_Labels})
    elif Host_Name == 'VictorChengPC':
        np.savez('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                 'Ave Photo Cntr Value with upper, mz_simed_out=%.3f, T=%.3eKData_Sel.npz' %
                 (mz_simed_out, Atom_Temperature), AveStd4All_Sel=AveStd4All_Sel, Data_Labels=Data_Labels)
        sio.savemat('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value with upper, mz_simed_out=%.3f, T=%.3eKData_Sel.mat' %
                    (mz_simed_out, Atom_Temperature), {'AveStd4All_Sel': AveStd4All_Sel, 'Data_Labels': Data_Labels})
    AveStd4All_Sel_flipped = np.delete(AveStd4All_Sel_flipped, for_del, 0)
    fig = plt.figure(figsize=(20, 15))
    # plt.plot(ky_flipped, Ave_Classify_All_Fit_flipped, 'r-', linewidth=5)
    plt.errorbar(AveStd4All_Sel_flipped[:, 8], AveStd4All_Sel_flipped[:, 0], AveStd4All_Sel_flipped[:, 1],
                 fmt='bo-', label='M', linewidth=5, markersize=20, mec='k', mfc='k')
    plt.errorbar(AveStd4All_Sel_flipped[:, 8], AveStd4All_Sel_flipped[:, 2], AveStd4All_Sel_flipped[:, 3],
                 fmt='mo-', label=r'$X_1$', linewidth=5, markersize=20, mec='k', mfc='k')
    plt.errorbar(AveStd4All_Sel_flipped[:, 8], AveStd4All_Sel_flipped[:, 4], AveStd4All_Sel_flipped[:, 5],
                 fmt='yo-', label=r'$X_2$', linewidth=5, markersize=20, mec='k', mfc='k')
    plt.errorbar(AveStd4All_Sel_flipped[:, 8], AveStd4All_Sel_flipped[:, 6], AveStd4All_Sel_flipped[:, 7],
                 fmt='go-', label=r'$\Gamma$', linewidth=5, markersize=20, mec='k', mfc='k')

    # plt.plot((np.min(ky_flipped), np.max(ky_flipped)), (0, 0))
    plt.title('Selected, SD')
    plt.xlabel(r'kz/$\pi$')
    plt.ylabel(r'$P(\Lambda_i)$')
    plt.xlim([np.min(ky_flipped), np.max(ky_flipped)])
    plt.ylim([-1, 1])

    # plt.plot(detuning_PWM_Weyl, 0, 'o', ms=10)
    # plt.title('Ave Photo Cntr Value with upper, Weyl:detuning=%.3f, T=%.3eK(SD)' %
    #           (detuning_PWM_Weyl, Atom_Temperature))
    plt.legend(loc='upper right')
    plt.tick_params(labelsize=30)
    plt.show()
    if Host_Name == 'calcenter2':
        fig.savefig('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value with upper, mz_simed_out=%.3f, T=%.3eK(SD)_Sel.png' %
                    (mz_simed_out, Atom_Temperature))
    elif Host_Name == 'VictorChengPC':
        fig.savefig('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value with upper, mz_simed_out=%.3f, T=%.3eK(SD)_Sel.png' %
                    (mz_simed_out, Atom_Temperature))

    LenSPs_Sel = np.array(LenSPs_Sel)
    LenSPs_flipped_Sel = np.concatenate((LenSPs_Sel, np.flip(LenSPs_Sel)))
    LenSPs_flipped_Sel = np.delete(LenSPs_flipped_Sel, for_del, 0)
    fig = plt.figure(figsize=(20, 15))
    # plt.plot(ky_flipped, Ave_Classify_All_Fit_flipped, 'r-', linewidth=5)
    plt.errorbar(AveStd4All_Sel_flipped[:, 8], AveStd4All_Sel_flipped[:, 0],
                 AveStd4All_Sel_flipped[:, 1] / np.sqrt(LenSPs_flipped_Sel),
                 fmt='bo-', label='M', linewidth=5, markersize=15, mec='k', mfc='k')
    plt.errorbar(AveStd4All_Sel_flipped[:, 8], AveStd4All_Sel_flipped[:, 2],
                 AveStd4All_Sel_flipped[:, 3] / np.sqrt(LenSPs_flipped_Sel),
                 fmt='mo-', label=r'$X_1$', linewidth=5, markersize=15, mec='k', mfc='k')
    plt.errorbar(AveStd4All_Sel_flipped[:, 8], AveStd4All_Sel_flipped[:, 4],
                 AveStd4All_Sel_flipped[:, 5] / np.sqrt(LenSPs_flipped_Sel),
                 fmt='yo-', label=r'$X_2$', linewidth=5, markersize=15, mec='k', mfc='k')
    plt.errorbar(AveStd4All_Sel_flipped[:, 8], AveStd4All_Sel_flipped[:, 6],
                 AveStd4All_Sel_flipped[:, 7] / np.sqrt(LenSPs_flipped_Sel),
                 fmt='go-', label=r'$\Gamma$', linewidth=5, markersize=15, mec='k', mfc='k')

    plt.plot((np.min(ky_flipped), np.max(ky_flipped)), (0, 0))
    plt.title('Selected, SE')
    plt.xlabel(r'kz/$\pi$')
    plt.ylabel(r'$P(\Lambda_i)$')
    plt.xlim([np.min(ky_flipped), np.max(ky_flipped)])
    plt.ylim([-1, 1])
    # plt.plot(detuning_PWM_Weyl, 0, 'o', ms=10)
    # plt.title('Ave Photo Cntr Value with upper, Weyl:detuning=%.3f, T=%.3eK(SE)' %
    #           (detuning_PWM_Weyl, Atom_Temperature))
    plt.legend(loc='upper right')
    plt.tick_params(labelsize=30)
    plt.show()
    if Host_Name == 'calcenter2':
        fig.savefig('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value with upper, mz_simed_out=%.3f, T=%.3eK(SE)_Sel.png' %
                    (mz_simed_out, Atom_Temperature))
    elif Host_Name == 'VictorChengPC':
        fig.savefig('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value with upper, mz_simed_out=%.3f, T=%.3eK(SE)_Sel.png' %
                    (mz_simed_out, Atom_Temperature))

    # Lower, sel
    # Pair_Data = sio.loadmat(Computer_Name + "Rb87Data\\3DSOCbandData\\190827\\"
    #                                         "detuning_ky_pair,mz_sim=%.3f.mat" % mz_simed_out)
    Detuning_ky_pair = Pair_Data['detuning_ky_pair']
    if mz_simed_out == -0.2:
        Detuning_ky_pair[0, 1] = 1
        Detuning_ky_pair[-1, 1] = 0
    for_del = np.asarray(Detuning_ky_pair[:, 1] == 100).nonzero()
    Detuning_ky_pair = np.delete(Detuning_ky_pair, for_del, 0)
    Detunings = Detuning_ky_pair[:, 0]
    kys = Detuning_ky_pair[:, 1]
    AveStd4All_Lower_Sel[:, 8] = kys
    kys_sort_arg = list(np.argsort(-kys))

    kys = kys[kys_sort_arg]
    Detunings = Detunings[kys_sort_arg]
    AveStd4All_Lower_Sel = AveStd4All_Lower_Sel[kys_sort_arg, :]

    ky_flipped = np.concatenate((kys, -np.flip(kys)))+1
    AveStd4All_Lower_Sel_flipped = \
        np.concatenate((np.flip(AveStd4All_Lower_Sel, axis=0), AveStd4All_Lower_Sel), axis=0)
    AveStd4All_Lower_Sel_flipped[:, 8] = ky_flipped
    for_del = np.asarray(AveStd4All_Lower_Sel_flipped[:, -1]).nonzero()
    if Host_Name == 'calcenter2':
        np.savez('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                 'Ave Photo Cntr Value no upper, mz_simed_out=%.3f, T=%.3eKData_Sel.npz' %
                 (mz_simed_out, Atom_Temperature), AveStd4All_Lower_Sel=AveStd4All_Lower_Sel, Data_Labels=Data_Labels)
        sio.savemat('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value no upper, mz_simed_out=%.3f, T=%.3eKData_Sel.mat' %
                    (mz_simed_out, Atom_Temperature),
                    {'AveStd4All_Lower_Sel': AveStd4All_Lower_Sel, 'Data_Labels': Data_Labels})
    elif Host_Name == 'VictorChengPC':
        np.savez('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                 'Ave Photo Cntr Value no upper, mz_simed_out=%.3f, T=%.3eKData_Sel.npz' %
                 (mz_simed_out, Atom_Temperature), AveStd4All_Lower_Sel=AveStd4All_Lower_Sel, Data_Labels=Data_Labels)
        sio.savemat('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value no upper, mz_simed_out=%.3f, T=%.3eKData_Sel.mat' %
                    (mz_simed_out, Atom_Temperature),
                    {'AveStd4All_Lower_Sel': AveStd4All_Lower_Sel, 'Data_Labels': Data_Labels})
    AveStd4All_Lower_Sel_flipped = np.delete(AveStd4All_Lower_Sel_flipped, for_del, 0)
    fig = plt.figure(figsize=(20, 15))
    plt.errorbar(AveStd4All_Lower_Sel_flipped[:, 8], AveStd4All_Lower_Sel_flipped[:, 0],
                 AveStd4All_Lower_Sel_flipped[:, 1],
                 fmt='bo-', label='M', linewidth=5, markersize=20, mec='k', mfc='k')
    plt.errorbar(AveStd4All_Lower_Sel_flipped[:, 8], AveStd4All_Lower_Sel_flipped[:, 2],
                 AveStd4All_Lower_Sel_flipped[:, 3],
                 fmt='mo-', label=r'$X_1$', linewidth=5, markersize=20, mec='k', mfc='k')
    plt.errorbar(AveStd4All_Lower_Sel_flipped[:, 8], AveStd4All_Lower_Sel_flipped[:, 4],
                 AveStd4All_Lower_Sel_flipped[:, 5],
                 fmt='yo-', label=r'$X_2$', linewidth=5, markersize=20, mec='k', mfc='k')
    plt.errorbar(AveStd4All_Lower_Sel_flipped[:, 8], AveStd4All_Lower_Sel_flipped[:, 6],
                 AveStd4All_Lower_Sel_flipped[:, 7],
                 fmt='go-', label=r'$\Gamma$', linewidth=5, markersize=20, mec='k', mfc='k')
    plt.plot((np.min(ky_flipped), np.max(ky_flipped)), (0, 0))
    plt.xlabel(r'kz/$\pi$')
    plt.ylabel(r'$P(\Lambda_i)$')
    plt.xlim([np.min(ky_flipped), np.max(ky_flipped)])
    plt.ylim([-1, 1])
    plt.legend(loc='upper right')
    plt.tick_params(labelsize=30)
    plt.title('Selected, SD, Lower')
    plt.show()

    fig = plt.figure(figsize=(20, 15))
    plt.errorbar(AveStd4All_Lower_Sel_flipped[:, 8], AveStd4All_Lower_Sel_flipped[:, 0],
                 AveStd4All_Lower_Sel_flipped[:, 1] / np.sqrt(LenSPs_flipped_Sel),
                 fmt='bo-', label='$M$', linewidth=5, markersize=15, mec='k', mfc='k')
    plt.errorbar(AveStd4All_Lower_Sel_flipped[:, 8], AveStd4All_Lower_Sel_flipped[:, 2],
                 AveStd4All_Lower_Sel_flipped[:, 3] / np.sqrt(LenSPs_flipped_Sel),
                 fmt='mo-', label=r'$X_1$', linewidth=5, markersize=15, mec='k', mfc='k')
    plt.errorbar(AveStd4All_Lower_Sel_flipped[:, 8], AveStd4All_Lower_Sel_flipped[:, 4],
                 AveStd4All_Lower_Sel_flipped[:, 5] / np.sqrt(LenSPs_flipped_Sel),
                 fmt='yo-', label=r'$X_2$', linewidth=5, markersize=15, mec='k', mfc='k')
    plt.errorbar(AveStd4All_Lower_Sel_flipped[:, 8], AveStd4All_Lower_Sel_flipped[:, 6],
                 AveStd4All_Lower_Sel_flipped[:, 7] / np.sqrt(LenSPs_flipped_Sel),
                 fmt='go-', label=r'$\Gamma$', linewidth=5, markersize=15, mec='k', mfc='k')
    plt.plot((np.min(ky_flipped), np.max(ky_flipped)), (0, 0))
    plt.xlabel(r'kz/$\pi$')
    plt.ylabel(r'$P(\Lambda_i)$')
    plt.xlim([np.min(ky_flipped), np.max(ky_flipped)])
    plt.ylim([-1, 1])
    plt.legend(loc='upper right')
    plt.tick_params(labelsize=30)
    plt.title('Selected, SE, Lower')
    plt.show()
    if Host_Name == 'calcenter2':
        fig.savefig('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value no upper, mz_simed_out=%.3f, T=%.3eK(SE)_Sel.png' %
                    (mz_simed_out, Atom_Temperature))
    elif Host_Name == 'VictorChengPC':
        fig.savefig('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                    'Ave Photo Cntr Value no upper, mz_simed_out=%.3f, T=%.3eK(SE)_Sel.png' %
                    (mz_simed_out, Atom_Temperature))

# %% Fitting with upper
# Detunings = np.linspace(-0.02, -0.6, 30)
Pair_Data = sio.loadmat(Computer_Name + "Rb87Data\\3DSOCbandData\\190827\\ky_detuning_pair,mz_sim=%.3f.mat" % mz_simed)
ky_Detuning_Pair = Pair_Data['ky_detuning_pair']
if Try == 4 and mz_simed == -0.25:
    ky_Detuning_Pair = np.flip(ky_Detuning_Pair)
Detunings = np.array(Detunings)

p_res, pcov = curve_fit(linear, Detunings, Ave_Classify_All, sigma=Std_Classify_All)
Ave_Classify_All_Fit = linear(Detunings, p_res[0], p_res[1])
kys = ky_Detuning_Pair[:, 1]
ky_flipped = np.concatenate((ky_Detuning_Pair[:, 1], -np.flip(ky_Detuning_Pair[:, 1])))
Ave_Classify_All_Fit_flipped = np.concatenate((Ave_Classify_All_Fit, np.flip(Ave_Classify_All_Fit)))
Ave_M_All_flipped = np.concatenate((Ave_M_All, np.flip(Ave_M_All)))
Std_M_All_flipped = np.concatenate((Std_M_All, np.flip(Std_M_All)))
Ave_X1_All_flipped = np.concatenate((Ave_X1_All, np.flip(Ave_X1_All)))
Std_X1_All_flipped = np.concatenate((Std_X1_All, np.flip(Std_X1_All)))
Ave_X2_All_flipped = np.concatenate((Ave_X2_All, np.flip(Ave_X2_All)))
Std_X2_All_flipped = np.concatenate((Std_X2_All, np.flip(Std_X2_All)))
Ave_Gamma_All_flipped = np.concatenate((Ave_Gamma_All, np.flip(Ave_Gamma_All)))
Std_Gamma_All_flipped = np.concatenate((Std_Gamma_All, np.flip(Std_Gamma_All)))

fig = plt.figure(figsize=(20, 15))
# plt.plot(ky_flipped, Ave_Classify_All_Fit_flipped, 'r-', linewidth=5)
plt.errorbar(ky_flipped, Ave_M_All_flipped, Std_M_All_flipped,
             fmt='bo-', label='M', linewidth=5, markersize=20, mec='k', mfc='k')
plt.errorbar(ky_flipped, Ave_X1_All_flipped, Std_X1_All_flipped,
             fmt='mo-', label=r'$X_1$', linewidth=5, markersize=20, mec='k', mfc='k')
plt.errorbar(ky_flipped, Ave_X2_All_flipped, Std_X2_All_flipped,
             fmt='yo-', label=r'$X_2$', linewidth=5, markersize=20, mec='k', mfc='k')
plt.errorbar(ky_flipped, Ave_Gamma_All_flipped, Std_Gamma_All_flipped,
             fmt='go-', label=r'$\Gamma$', linewidth=5, markersize=20, mec='k', mfc='k')

plt.plot((np.min(ky_flipped), np.max(ky_flipped)), (0, 0))
# plt.title('Average Photo Center Value')
plt.xlabel(r'kz/$\pi$')
plt.ylabel(r'$P(\Lambda_i)$')
plt.xlim([np.min(ky_flipped), np.max(ky_flipped)])
plt.ylim([-0.5, 1])
detuning_PWM_Weyl = -p_res[1] / p_res[0]
# plt.plot(detuning_PWM_Weyl, 0, 'o', ms=10)
# plt.title('Ave Photo Cntr Value with upper, Weyl:detuning=%.3f, T=%.3eK(SD)' % (detuning_PWM_Weyl, Atom_Temperature))
plt.legend(loc='upper right')
plt.tick_params(labelsize=30)
plt.show()
os.system('mkdir ' + 'X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\')
if Host_Name == 'calcenter2':
    fig.savefig('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                'Ave Photo Cntr Value with upper, Weyl-detuning=%.3f, T=%.3eK(SD).png' %
                (detuning_PWM_Weyl, Atom_Temperature))
elif Host_Name == 'VictorChengPC':
    fig.savefig('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                'Ave Photo Cntr Value with upper, Weyl-detuning=%.3f, T=%.3eK(SD).png' %
                (detuning_PWM_Weyl, Atom_Temperature))

LenSPs = np.array(LenSPs)
LenSPs_flipped = np.concatenate((LenSPs, np.flip(LenSPs)))
fig = plt.figure(figsize=(20, 15))
# plt.plot(ky_flipped, Ave_Classify_All_Fit_flipped, 'r-', linewidth=5)
plt.errorbar(ky_flipped, Ave_M_All_flipped, Std_M_All_flipped / np.sqrt(LenSPs_flipped),
             fmt='bo-', label='M', linewidth=5, markersize=15, mec='k', mfc='k')
plt.errorbar(ky_flipped, Ave_X1_All_flipped, Std_X1_All_flipped / np.sqrt(LenSPs_flipped),
             fmt='mo-', label=r'$X_1$', linewidth=5, markersize=15, mec='k', mfc='k')
plt.errorbar(ky_flipped, Ave_X2_All_flipped, Std_X2_All_flipped / np.sqrt(LenSPs_flipped),
             fmt='yo-', label=r'$X_2$', linewidth=5, markersize=15, mec='k', mfc='k')
plt.errorbar(ky_flipped, Ave_Gamma_All_flipped, Std_Gamma_All_flipped / np.sqrt(LenSPs_flipped),
             fmt='go-', label=r'$\Gamma$', linewidth=5, markersize=15, mec='k', mfc='k')

plt.plot((np.min(ky_flipped), np.max(ky_flipped)), (0, 0))
# plt.title('Average Photo Center Value')
plt.xlabel(r'kz/$\pi$')
plt.ylabel(r'$P(\Lambda_i)$')
plt.xlim([np.min(ky_flipped), np.max(ky_flipped)])
plt.ylim([-0.5, 1])
detuning_PWM_Weyl = -p_res[1] / p_res[0]
# plt.plot(detuning_PWM_Weyl, 0, 'o', ms=10)
# plt.title('Ave Photo Cntr Value with upper, Weyl:detuning=%.3f, T=%.3eK(SE)' % (detuning_PWM_Weyl, Atom_Temperature))
plt.legend(loc='upper right')
plt.tick_params(labelsize=30)
plt.show()
if Host_Name == 'calcenter2':
    fig.savefig('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                'Ave Photo Cntr Value with upper, Weyl-detuning=%.3f, T=%.3eK(SE).png' %
                (detuning_PWM_Weyl, Atom_Temperature))
elif Host_Name == 'VictorChengPC':
    fig.savefig('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                'Ave Photo Cntr Value with upper, Weyl-detuning=%.3f, T=%.3eK(SE).png' %
                (detuning_PWM_Weyl, Atom_Temperature))

# %% Selected with upper
# %% Fitting no upper
# fig = plt.figure(figsize=(30, 15))
# plt.errorbar(Detunings, Ave_Lower_Classify_All, Std_Lower_Classify_All, fmt='*')
# plt.plot((np.min(Detunings), np.max(Detunings)), (0, 0))
# # plt.title('Average Photo Center Value')
# plt.xlabel('detuning/Er')
# plt.ylabel('Spin Polarization')
# plt.xlim([np.min(Detunings), np.max(Detunings)])
# plt.ylim([-1, 1])
# plt.title('Ave Photo Cntr Value without upper, T=%.3eK(SD)' % Atom_Temperature)
# plt.show()
AveStd4All_Lower = np.zeros((len(Ave_M_Lower_All), 9))
AveStd4All_Lower[:, 0] = Ave_M_Lower_All
AveStd4All_Lower[:, 1] = Std_M_Lower_All
AveStd4All_Lower[:, 2] = Ave_X1_Lower_All
AveStd4All_Lower[:, 3] = Std_X1_Lower_All
AveStd4All_Lower[:, 4] = Ave_X2_Lower_All
AveStd4All_Lower[:, 5] = Std_X2_Lower_All
AveStd4All_Lower[:, 6] = Ave_Gamma_Lower_All
AveStd4All_Lower[:, 7] = Std_Gamma_Lower_All
AveStd4All_Lower[:, 8] = kys
if Host_Name == 'calcenter2':
    np.savez('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
             'Ave Photo Cntr Value no upper, mz_simed_out=%.3f, T=%.3eKData.npz' %
             (mz_simed_out, Atom_Temperature), AveStd4All_Lower=AveStd4All_Lower, Data_Labels=Data_Labels)
    sio.savemat('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                'Ave Photo Cntr Value no upper, mz_simed_out=%.3f, T=%.3eKData.mat' %
                (mz_simed_out, Atom_Temperature),
                {'AveStd4All_Lower': AveStd4All_Lower, 'Data_Labels': Data_Labels})
elif Host_Name == 'VictorChengPC':
    np.savez('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
             'Ave Photo Cntr Value no upper, mz_simed_out=%.3f, T=%.3eKData.npz' %
             (mz_simed_out, Atom_Temperature), AveStd4All_Lower=AveStd4All_Lower, Data_Labels=Data_Labels)
    sio.savemat('X:\\work\\mail.ustc.edu.cn\\3DSOC_spin_texture\\figures\\Selected_Pic\\'
                'Ave Photo Cntr Value no upper, mz_simed_out=%.3f, T=%.3eKData.mat' %
                (mz_simed_out, Atom_Temperature),
                {'AveStd4All_Lower': AveStd4All_Lower, 'Data_Labels': Data_Labels})
Ave_M_Lower_All_flipped = np.concatenate((Ave_M_Lower_All, np.flip(Ave_M_Lower_All)))
Std_M_Lower_All_flipped = np.concatenate((Std_M_Lower_All, np.flip(Std_M_Lower_All)))
Ave_X1_Lower_All_flipped = np.concatenate((Ave_X1_Lower_All, np.flip(Ave_X1_Lower_All)))
Std_X1_Lower_All_flipped = np.concatenate((Std_X1_Lower_All, np.flip(Std_X1_Lower_All)))
Ave_X2_Lower_All_flipped = np.concatenate((Ave_X2_Lower_All, np.flip(Ave_X2_Lower_All)))
Std_X2_Lower_All_flipped = np.concatenate((Std_X2_Lower_All, np.flip(Std_X2_Lower_All)))
Ave_Gamma_Lower_All_flipped = np.concatenate((Ave_Gamma_Lower_All, np.flip(Ave_Gamma_Lower_All)))
Std_Gamma_Lower_All_flipped = np.concatenate((Std_Gamma_Lower_All, np.flip(Std_Gamma_Lower_All)))

fig = plt.figure(figsize=(20, 15))
plt.errorbar(ky_flipped, Ave_M_Lower_All_flipped, Std_M_Lower_All_flipped,
             fmt='bo-', label='M', linewidth=5, markersize=20, mec='k', mfc='k')
plt.errorbar(ky_flipped, Ave_X1_Lower_All_flipped, Std_X1_Lower_All_flipped,
             fmt='mo-', label=r'$X_1$', linewidth=5, markersize=20, mec='k', mfc='k')
plt.errorbar(ky_flipped, Ave_X2_Lower_All_flipped, Std_X2_Lower_All_flipped,
             fmt='yo-', label=r'$X_2$', linewidth=5, markersize=20, mec='k', mfc='k')
plt.errorbar(ky_flipped, Ave_Gamma_Lower_All_flipped, Std_Gamma_Lower_All_flipped,
             fmt='go-', label=r'$\Gamma$', linewidth=5, markersize=20, mec='k', mfc='k')
plt.plot((np.min(ky_flipped), np.max(ky_flipped)), (0, 0))
plt.xlabel(r'kz/$\pi$')
plt.ylabel(r'$P(\Lambda_i)$')
plt.xlim([np.min(ky_flipped), np.max(ky_flipped)])
plt.ylim([-1, 1])
plt.legend(loc='upper right')
plt.tick_params(labelsize=30)
plt.show()

fig = plt.figure(figsize=(20, 15))
plt.errorbar(ky_flipped, Ave_M_Lower_All_flipped, Std_M_Lower_All_flipped / np.sqrt(LenSPs_flipped),
             fmt='bo-', label='M', linewidth=5, markersize=15, mec='k', mfc='k')
plt.errorbar(ky_flipped, Ave_X1_Lower_All_flipped, Std_X1_Lower_All_flipped / np.sqrt(LenSPs_flipped),
             fmt='mo-', label=r'$X_1$', linewidth=5, markersize=15, mec='k', mfc='k')
plt.errorbar(ky_flipped, Ave_X2_Lower_All_flipped, Std_X2_Lower_All_flipped / np.sqrt(LenSPs_flipped),
             fmt='yo-', label=r'$X_2$', linewidth=5, markersize=15, mec='k', mfc='k')
plt.errorbar(ky_flipped, Ave_Gamma_Lower_All_flipped, Std_Gamma_Lower_All_flipped / np.sqrt(LenSPs_flipped),
             fmt='go-', label=r'$\Gamma$', linewidth=5, markersize=15, mec='k', mfc='k')
plt.plot((np.min(ky_flipped), np.max(ky_flipped)), (0, 0))
plt.xlabel(r'kz/$\pi$')
plt.ylabel(r'$P(\Lambda_i)$')
plt.xlim([np.min(ky_flipped), np.max(ky_flipped)])
plt.ylim([-1, 1])
plt.legend(loc='upper right')
plt.tick_params(labelsize=30)
plt.show()
# fig.savefig(r'D:\Rb87figures\3DSOC_spin_texture\figures\Selected_Pic\CenterValueMeta%dto%d.png' % (
#     Meta_Row_List[0], Meta_Row_List[1]))
# fig.savefig(r'D:\Rb87figures\3DSOC_spin_texture\figures\Selected_Pic\CornerValueMeta33to76.png')

# %% Calculate selected
'''Data_Table = xlrd.open_workbook(Computer_Name + 'Rb87files\\3DSOC_spin_texture\\spin_texture.xlsx').sheets()[0]
Meta_Table = xlrd.open_workbook(Computer_Name + 'Rb87files\\3DSOC_spin_texture\\spin_texture.xlsx').sheets()[1]
Fitted = Ave_Center_All_Fit
Selected_File = open(Computer_Name + "Rb87Data\\3DSOC_spin_texture\\Selected\\"
                                     "SelectedMeta_Row=%d,Detuning=%.3f.txt" % (Meta_Row, Detuning), 'w+')
Ave_Center_All_Selected = []
Detunings_Selected = []
for Meta_Row_Index in range(len(Meta_Row_List)):
    Meta_Row = Meta_Row_List[Meta_Row_Index]
    Meta = Meta_Table.row_values(Meta_Row - 1)
    Try = int(Meta[1])
    Detuning = Meta[2]
    Mat_Data_command = Meta[3]
    exec(Mat_Data_command)
    Bad_Dict_Command = Meta[4]
    exec(Bad_Dict_Command)
    SP_Selected = {}
    SPs = {}
    for Line_Index in range(len(Mat_Data)):
        Line = Mat_Data[Line_Index]
        Input_Matrix = Data_Table.row_values(Line - 1)
        YearMonthDay = str(int(Input_Matrix[0]))
        Data_Folder_Name = Computer_Name + 'Rb87Data\\3DSOC_spin_texture\\' + \
            Fit_Prefix[Fit_BEC_Switch] + '\\' + YearMonthDay + '\\All_Spin_Polarization\\'
        Ini_Photo = int(Input_Matrix[1])
        Fin_Photo = int(Input_Matrix[2])
        if Line in Bad_Dict.keys():
            Bad_List = Bad_Dict[Line]
        else:
            Bad_List = []
        for Photo_Count in range(Fin_Photo - Ini_Photo + 1):
            Photo_Ind = Photo_Count + Ini_Photo
            if Photo_Ind in Bad_List:
                continue
            Data = np.load(Data_Folder_Name + 'Cam_22860_photo%d_Spin_Polarization_FBZ.npz' % Photo_Ind)
            SP = Data['SP']
            if Detuning < 0:
                Classify = m_data
            else:
                Classify = gamma_data
            Indicator = Classify(SP, Size=Classify_Size, Find_Center=1)[0]
            SPs[(Line, Photo_Ind)] = SP
            if np.abs(Indicator - Fitted[Meta_Row_Index]) < 0.05:
                SP_Selected[(Line, Photo_Ind)] = SP

    print("Row %d Selected %d out of %d photos, Detuning=%.3f." % (Meta_Row, len(SP_Selected), len(SPs), Detuning))
    Selected_File.write("Row %d Selected %d out of %d photos, Detuning=%.3f.\n" %
                        (Meta_Row, len(SP_Selected), len(SPs), Detuning))
    SP_Ave_Selected = np.average(np.array(list(SP_Selected.values())), axis=0)
    SP_Std_Selected = np.std(np.array(list(SP_Selected.values())), axis=0)
    if np.isnan(SP_Ave_Selected).any():
        pass
    else:
        Detunings_Selected.append(Detuning)
        Ave_Center_All_Selected.append(Classify(SP_Ave_Selected, Size=Classify_Size, Find_Center=1)[0])

    np.savez(Computer_Name + "Rb87Data\\3DSOC_spin_texture\\Selected\\"
             "ClassifyMeta_Row=%d,Detuning=%.3f_Selected.npz" % (Meta_Row, Detuning),
             Indicators=Indicators, Mat_Data=Mat_Data, Bad_Dict=Bad_Dict, SP_Ave=SP_Ave_Selected)
    sio.savemat(Computer_Name + "Rb87Data\\3DSOC_spin_texture\\Selected\\"
                "ClassifyMeta_Row=%d,Detuning=%.3f_Selected.mat" % (Meta_Row, Detuning),
                {'Indicators': Indicators, 'Mat_Data': Mat_Data, 'Bad_Dict': Bad_Dict, 'SP_Ave': SP_Ave_Selected})

Selected_File.close()


# %%
Detunings_Selected = np.array(Detunings_Selected)
p_res_Selected, pcov = curve_fit(linear, Detunings_Selected, Ave_Center_All_Selected)#, sigma=SP_Std_Selected)
Ave_Center_All_Selected_Fit = linear(Detunings_Selected, p_res_Selected[0], p_res_Selected[1])

fig = plt.figure(1001)
plt.plot(Detunings_Selected, Ave_Center_All_Selected, '*', Detunings_Selected, Ave_Center_All_Selected_Fit, '-')
plt.plot((np.min(Detunings_Selected), np.max(Detunings_Selected)), (0, 0))

plt.xlabel('detuning/Er')
plt.ylabel('Spin Polarization')
plt.xlim([np.min(Detunings_Selected), np.max(Detunings_Selected)])
plt.ylim([np.min((Ave_Center_All_Selected, Ave_Center_All_Selected_Fit)),
          np.max((Ave_Center_All_Selected, Ave_Center_All_Selected_Fit))])

detuning_PWM_Weyl_Selected = -p_res_Selected[1] / p_res_Selected[0]
plt.plot(detuning_PWM_Weyl_Selected, 0, 'o', ms=10)
plt.title('Average Photo Center Value Selected, Weyl=%.3f' % detuning_PWM_Weyl_Selected)
plt.show(block=False)
plt.pause(0.5)
plt.close(1001)
fig.savefig('D:\\Rb87figures\\3DSOC_spin_texture\\figures\\Selected_Pic\\CenterValueMeta%dto%d_Selected.png' % (
    Meta_Row_List[0], Meta_Row_List[1]))

np.savez(Computer_Name + "Rb87Data\\3DSOC_spin_texture\\Selected\\"
         "ClassifyMeta_RowEnd=%d,Classifys.npz" % Meta_Row,
         Detunings_Selected=Detunings_Selected, Ave_Center_All_Selected=Ave_Center_All_Selected,
         Ave_Center_All_Selected_Fit=Ave_Center_All_Selected_Fit,
         detuning_PWM_Weyl_Selected=detuning_PWM_Weyl_Selected)

sio.savemat(Computer_Name + "Rb87Data\\3DSOC_spin_texture\\Selected\\"
            "ClassifyMeta_RowEnd=%d,Classifys" % Meta_Row,
            {'Detunings_Selected': Detunings_Selected, 'Ave_Center_All_Selected': Ave_Center_All_Selected,
             'Ave_Center_All_Selected_Fit': Ave_Center_All_Selected_Fit,
             'detuning_PWM_Weyl_Selected': detuning_PWM_Weyl_Selected})'''

#%%
np.savez("C:\\Users\\cheng\\Desktop\\pht\\data.npz",
         kz=ky_flipped,
         M=Ave_M_All_flipped, M_Std=Std_M_All_flipped,
         X1=Ave_X1_All_flipped, X1_Std=Std_X1_All_flipped,
         X2=Ave_X2_All_flipped, X2_Std=Std_X2_All_flipped,
         Gamma=Ave_Gamma_All_flipped, Gamma_Std=Std_Gamma_All_flipped,
         M_Lower=Ave_M_Lower_All_flipped, M_Lower_Std=Std_M_Lower_All_flipped,
         X1_Lower=Ave_X1_Lower_All_flipped, X1_Lower_Std=Std_X1_Lower_All_flipped,
         X2_Lower=Ave_X2_Lower_All_flipped, X2_Lower_Std=Std_X2_Lower_All_flipped,
         Gamma_Lower=Ave_Gamma_Lower_All_flipped, Gamma_Lower_Std=Std_Gamma_Lower_All_flipped)
Data = np.load("C:\\Users\\cheng\\Desktop\\pht\\data.npz")
Data_dict = dict(Data)
sio.savemat("C:\\Users\\cheng\\Desktop\\pht\\data.mat", Data_dict)
time_end = time.time()
print('Time spent: ', time_end - time_start, 's')
time_end_struct = time.localtime(time_end)
print('Start time:', time_start_struct)
print('End time:', time_end_struct)
